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1  Introduction 

Consider  a  patch  of  dye  immersed  in  a  flowing  fluid,  confined  within  a  spatial  domain  0. 
The  concentration  of  the  dye  9(x,t )  can  be  modelled  by  the  advection-diffusion  equation 

c)Q 

—  +u-V9  =  kV29,  (1) 

where  n  is  a  coefficient  of  diffusivity  and  u(x,t )  is  a  prescribed  velocity  field.  By  using 
equation  (1),  it  is  assumed  that  the  concentration  is  advected  precisely  with  the  flow,  and 
is  also  subject  to  molecular  diffusion  by  which  the  concentration  is  allowed  to  gradually 
spread  to  adjacent  fluid  elements.  The  velocity  field  in  equation  (1)  is  assumed  to  have 
been  determined  a  priori  and  satisfies  the  incompressibility  condition 


V-w  =  0.  (2) 

We  also  assume  that  the  boundary  conditions 

n  ■  u  =  0,  (3) 

n  ■  V0  =  0,  (4) 

apply  on  the  boundary  of  the  domain  <9R,  where  n  is  the  unit  outward  normal  to  the 
boundary.  These  conditions  ensure  that  no  fluid  penetrates  the  boundary  of  the  domain 
and  that  no  concentration  is  allowed  to  diffuse  into  or  out  of  the  domain,  respectively.  (We 
can  also  use  spatially-periodic  boundary  conditions,  which  ensure  that  there  is  no  net  flux 
in  or  out  of  the  domain.) 

We  proceed  to  demonstrate  some  general  properties  of  the  system  defined  by  equations 
(1)— (4)  by  deriving  evolution  equations  for  the  spatial  mean  and  variance  of  9.  To  do  this 
we  begin  by  multiplying  equation  (1)  by  m9m~ 1  and  reordering  the  derivatives  to  give 

^ :(9m )  +  V  •  ( 9mu )  =  urn  [V  •  (9m~LV9)  -  (m  -  l)9m~2\V  9\2]  .  (5) 

Now  we  take  the  volume  integral  over  the  domain,  defined  by 

(4)  =  [  4>  dV,  (6) 
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where  cf>(x,t)  is  any  function,  which  gives 

Jt{ 9m )  =  —m(m  -  l)K(9m~2\V9\2), 

on  applying  the  divergence  theorem  to  the  second  term  on  the  left-hand  side  of  equation 
(5)  and  the  first  term  on  the  right-hand  side,  and  applying  boundary  conditions  (3)-(4)  to 
make  the  resulting  boundary  terms  vanish.  Setting  m  =  1,2  gives  the  two  equations 

>  =  <8> 

|(02)  =  -2K(|V0|2),  (9) 

respectively.  The  first  equation  above  represents  conservation  of  total  concentration,  which 
is  expected  given  that  no  concentration  enters  the  domain  as  ensured  by  boundary  conditions 
(3)-(4).  The  second  equation  is  an  evolution  equation  for  the  second  moment  of  9  in  terms 
of  its  gradient.  This  equation  is  not  closed,  i.e.  it  does  not  independently  determine  the  time 
evolution  of  (92),  because  its  right-hand  side  depends  on  |V0|,  whose  time  evolution  must 
be  determined  from  a  solution  of  the  advection-diffusion  equation  (1).  However,  noting  that 
the  right-hand  side  of  equation  (9)  is  negative-definite  if  9  is  not  constant  in  space,  we  can 
deduce  that  the  second  moment  of  the  concentration  must  decrease  in  time.  Thus,  at  large 
times  the  second  moment  converges  toward  a  constant  value  in  which  the  right-hand  side  of 
(9)  vanishes.  This  occurs  if  and  only  if  the  gradient  of  the  concentration  is  zero  everywhere. 
The  system  is  driven  toward  a  fully  homogenized  state  in  which  the  second  moment  takes 
the  value 

(92)  =  (10) 

Note  that  although  (92)  never  quite  reaches  this  value,  we  say  the  system  is  mixed  once  it 
has  fallen  below  a  certain  prescribed  threshold. 

As  mentioned  previously,  equation  (9)  does  not  independently  determine  the  rate  at 
which  (92)  decreases  in  time.  In  fact,  the  rate  at  which  (92)  decreases  usually  depends 
strongly  on  the  prescribed  velocity  field,  which  does  not  appear  in  equation  (9).  Importantly, 
the  introduction  of  a  velocity  field,  which  we  refer  to  as  stirring ,  usually  has  the  effect  of 
increasing  gradients  in  the  concentration  field  |V#|  through  advection  and  thus  has  the 
potential  to  significantly  accelerate  the  onset  of  mixing. 


filaments 


Figure  1:  Typical  mixing  scenario  in  which  a  patch  of  dye  is  advected  into  filaments. 
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Figure  2:  Decay  of  the  second  moment  of  the  concentration  field  under  the  action  of  stirring. 

A  typical  stirring  scenario  is  illustrated  in  figure  1,  in  which  an  initially  isolated  patch  of 
dye  is  advected  into  thin  striations  that  exhibit  a  much  greater  concentration  gradient  than 
the  initial  patch.  The  evolution  of  the  second  moment,  corresponding  to  this  typical  mixing 
scenario,  is  shown  in  figure  2.  During  the  initial  stirring  phase,  before  thin  filaments  are 
formed,  the  concentration  varies  over  relatively  large  length  scales  and  the  rate  of  mixing 
will  be  slow.  After  a  short  transient,  by  which  point  the  filaments  have  become  sufficiently 
thin,  the  rate  of  diffusion  is  large  enough  to  balance  the  rate  of  thinning  of  the  filaments 
and  the  rate  of  mixing  d(82)/dt  is  no  longer  characterized  by  n,  but  is  instead  characterized 
by  the  rate  at  which  the  velocity  field  is  straining  the  fluid  into  filaments.  Once  this  occurs 
the  characteristic  width  of  the  filaments  remains  constant,  and  the  rate  of  mixing  reaches 
a  maximum  which  persists  for  all  time.  The  independence  of  the  rate  of  diffusion  on  the 
parameter  n  during  the  mixing  phase  implies  the  scaling 

|V0|  ~  ft- 1/2 ,  (11) 

which  can  be  used  as  an  indication  of  good  mixing  in  practical  applications. 

The  coefficient  of  diffusion  k  is  small  in  most  applications.  Hence,  the  characteristic 
timescale  over  which  diffusion  takes  place  in  the  absence  of  an  advective  flow  (u  =  0) 
is  typically  very  long.  For  example,  in  the  context  of  heat  diffusing  through  air,  k  ~ 
2.2  x  10~5m2/s,  so  in  a  room  with  sides  of  length  L  =  10m,  a  characteristic  diffusion 
timescale  is  L2 / k  ~  53  days.  Thus,  if  one  is  aiming  to  mix  the  concentration  field,  relying 
on  diffusion  in  the  absence  of  advection  is  typically  ineffective.  In  the  following  section  we 
present  an  example  that  can  be  solved  analytically  to  demonstrate  the  development  of  the 
mixing  phase  discussed  above. 
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The  filament  solution 


Consider  the  prescribed  velocity  field 

u(x)  =  (Ax, -A  y), 


(12) 


where  A  is  a  constant  parameter,  which  drives  fluid  towards  the  x-axis  in  a  straining  motion 
with  a  hyperbolic  point  at  the  origin.  In  this  case  the  advection-diffusion  equation  takes 
the  form 

80  .89  80  2 

at  +  Xxai  ~  x%  =  kV  (13) 

This  equation  can  be  solved  exactly,  but  for  simplicity  here  we  consider  only  x-independent 
solutions  of  the  form 

0(y,  t)  =  e~xtf(y).  (14) 

Substituting  this  into  equation  (13)  gives  the  ordinary  differential  equation 

-A/  -  A yf  =  nf" .  (15) 


If  we  further  assume  that  the  concentration  vanishes  in  the  far  field,  /  — >  0  as  y  — >  ±oo,  we 
can  solve  equation  (15)  and  derive  the  solution 


e  =  e-xte-y2/212, 


(16) 


where 

i  =  (17) 

is  called  the  Batchelor  length. 

The  solution  (16)  represents  a  filament  that  is  an  attractor  for  typical  initial  conditions 
of  the  system.  For  any  reasonable  initial  concentration  field,  the  compression  due  to  the 
straining  velocity  will  attempt  to  confine  the  concentration  into  an  ever  thinner  region  about 
the  x-axis  until  its  gradient  becomes  sufficiently  large  for  diffusion  to  balance  the  rate  of 
compression.  The  width  of  the  resulting  filament  is  characterized  by  the  Batchelor  length 
(17),  which  contains  the  necessary  scaling  for  (11)  to  apply.  Furthermore,  we  see  from  the 
solution  (16)  that  the  rate  of  mixing  A  during  the  mixing  phase  is  indeed  independent  of  k. 
Note  that  equation  (8)  must  be  taken  with  a  grain  of  salt,  because  the  concentration  does 
not  vanish  as  |x|  — >  oo.  This  is  a  consequence  of  the  neglect  of  the  x  dependence. 


2  Effective  Diffusivity 

Recall  that  filaments  undergoing  advection  are  strained  along  one  direction  at  a  rate  A  and 
are  compressed  along  the  other  to  the  Batchelor  length  l  =  \J k/ A,  where  n  is  the  molecular 
diffusivity  and  A  is  the  local  rate  of  strain.  In  the  previous  section  it  was  shown  that  the 
variance  evolves  according  to 

(O'2)  ~  e~xt , 

where  A  is  a  local  rate  of  strain. 

From  simple  arguments,  molecular  diffusivity  alone  is  an  exceedingly  weak  method  for 
mixing  scalar  fields.  We  may  then  ask  how  stirring  might  improve  tracer  mixing  and,  more 
importantly,  how  we  can  quantify  the  effects  of  stirring  by  invoking  an  effective  diffusivity. 
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2.1  A  random  walk  model 

The  advection-diffusion  equation  for  a  passive  scalar  8  =  9(x,  t )  can  be  written  as 


—  +  u  ■  V9  =  k\729,  (18) 

ot 

where  u  is  a  given  velocity  field  and  k  is  a  constant  of  diffusivity.  Even  in  the  absence  of 
stirring  and  the  hlamentation  process,  this  diffusivity  can  act  to  mix  the  tracer  field  8  by 
monotonously  decreasing  the  variance  from  equation  (9). 

The  long  timescales  associated  with  diffusion  alone  indicates  that  the  physical  act  of 
stirring,  or  the  advective  term  u  ■  V 8  of  (18),  is  necessary  to  efficiently  mix  a  tracer  held. 
However,  the  advective  term  of  (18)  is  more  complex  than  the  diffusivity.  Ideally,  one 
could  write  the  effects  of  both  the  mixing  of  the  tracer  held  by  molecular  diffusivity  and 
by  advective  stirring  as  an  effective  diffusivity,  such  that  the  equations  of  motion  can  be 
approximated  by 

BO 

—  =  K\72e  -u-V8k  KeV28.  (19) 

We  consider  the  displacement,  x,  of  an  individual  particle  stirred  by  a  random  how  held. 
We  would  like  to  compare  the  displacemnt  of  the  particle  to  a  true  random  walk, 

Xn  —  Xn—  1  T  £n  (20) 

where  xn  is  position  of  the  particle  at  time  n  and  where  £ n  is  a  random  kick.  Here  we 
assume  that  all  £n  are  independent  and  identically  distributed  random  variables. 

At  the  nth  step,  the  particle  position  xn  satishes 

n 

xn  =  X0  + ^2  tk 

k= 1 

where  xq  is  the  initial  position.  If  the  initial  position  of  the  particle  is  xq  =  0  and  the  mean 
individual  kick  (£n)  =  0,  then  the  mean  position  of  the  particle  at  time  n  is 

n 

{ Xn )  =  ^(Zk)  =  0 
k=  1 

and  the  action  of  the  random  variations  in  the  how  held  given  by  £  do  not  change  the  mean 
position  of  the  particle. 

Consider  instead  the  mean-squared  displacement  of  the  particle: 

n 

{ xl )  =  ^2=  n(f)  =  n(j2’  (21) 

k= 1 

where  cr2  =  (£2).  The  mean-squared  displacement  of  the  particle  then  grows  linearly  with 
the  number  of  steps,  n,  at  a  rate  determined  by  the  standard  deviation  of  the  random 
forcing  held. 
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For  a  process  with  an  average  length  of  time  between  kicks,  T,  the  number  of  kicks  can 
be  converted  into  a  time  by  letting  t  =  nT.  The  squared  displacement  in  (21)  then  grows 
linearly  in  time  as 

(: x2(t ))  =  2  Kt  (22) 

where  k  =  a2  /  2 T  is  a  diffusivity,  or  a  measure  of  how  quickly  the  squared  displacement  of 
the  particle  grows  in  time. 

The  above  calculations  are  only  in  one  dimension.  In  d  dimensions,  the  mean-squared 
distance  from  the  origin  can  be  written  as  a  sum  of  the  mean-squared  distance  (xl  )  along 
each  dimension  k  at  timestep  n, 

d 

(^^xt,n)  =  (X1  ,n)  +  (x2,n)  +  •••  +  (xd,n) 
k= 1 

which,  assuming  the  random  forcing  is  isotropic  such  that  (x\  n)  =  (x2),  can  also  be  written 
as 


d(x2)  =  dna2  =  2  dnt  (23) 

where  k  =  a2 /2T  is  the  same  as  in  the  one  dimensional  case. 

The  diffusivity,  k,  can  be  thought  of  as  an  effective  diffusivity  as  it  arises  from  the 
physical  stirring  generated  by  the  stochastic  process  £n.  The  effective  diffusivity  is  again 
given  by 


(24) 


where  a2  is  the  standard  deviation  of  the  particle  displacement  due  to  the  random  kicks  £n 
and  T  is  the  average  interval  between  kicks. 


^eff 


a 

2 T 


2.2  Particle  concentration  and  diffusivity 

The  random  walk  description  of  a  diffusive  process  above  applies  only  to  a  single  particle 
with  position  given  by  xn.  We  wish  to  consider  the  evolution  of  a  cloud  of  particles.  Let 
6{x.  t )  be  the  density  of  particles  within  a  given  area  (or  volume  in  3D);  then  9  satisfies  the 
diffusive  equation 

=  Ke\720  (25) 

dt  y  ’ 

if  each  point  evolves  independently  according  to  a  random  walk  xn+\  =  xn  +  £n. 

The  idea  of  a  cloud  of  particles  only  works  if  the  spatial  and  temporal  resolution  of  the 
model  are  large  enough  such  that  the  motion  of  individual  particles  are  decorrelated.  If 
we  are  zoomed  in  too  closely,  such  that  the  grid  scale  is  less  than  the  standard  deviation 
of  the  displacements  ( L  <  a),  then  we  see  only  individual  particles  rather  than  a  cloud  of 
particles  with  a  given  concentration,  9.  If  we  choose  sampling  timescales  that  are  too  short, 
such  that  the  timescale  of  the  flow  is  less  than  the  interval  between  kicks,  T,  then  the  kicks 
become  correlated  in  time. 
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2.3  Numerical  examples  for  sine  flow 

Consider,  for  example,  the  advection  of  a  cloud  of  particles  by  the  famous  sine  flow  in 
a  doubly-periodic  domain  [2].  The  velocity  field  of  the  sine  flow  alternates  between  a 
horizontally-aligned  shear  flow  and  vertically-aligned  shear  flow  given  by 

uh  =  (U  sin  (27 xky/L)  ,0)  for  0  <  t  <  ^r; 
uy  =  (0  ,  U  sin  (2irkx/L))  for  \t  <t  <  t, 

where  here  U  is  a  constant  velocity,  k  is  the  wavenumber  of  the  shear  flow  field,  L  is  the 
size  of  the  biperiodic  domain,  and  r  is  the  period  of  the  flow. 

We  then  obtain  a  map  of  the  position  of  the  particle  by  solving 


X  =  U,  X  =  Xq 


exactly  as  a  two-step  process: 

Step  1  :  x(r/2 )  =  xq  +  U\  sin  ^  2n^yo  j 

y(r/2)  =  y0 


Step  2  :  x(r)  =  x{t / 2) 

y(T)  =  y[r/ 2)  +  U\  sin  (27r^r/2))  . 

At  the  end  of  an  interval  of  length  r,  the  new  position  (a of  a  particle  that  started  at 
(x,  y)  is  then 


x  =  x  +  T  sin 


27 rky\ 

~T~  ) 


y  =  y  +  T  sin 


2irkx'  \ 

) 


(26) 


where  T  =  Ut/2.  For  a  constant  velocity,  U,  T  is  proportional  to  r  and  the  intensity  of 
stirring  can  be  manipulated  by  varying  T.  Note  that  in  the  second  step  of  (26),  y' ,  is  a 
function  of  x'  and  not  x.  This  is  necessary  to  conserve  area,  as  can  be  verified  by  computing 
the  Jacobian  of  the  map  (26). 


2.4  Diffusivity  and  the  limit  of  small  T 

We  can  vary  the  period  of  the  sine  flow  given  in  (26)  by  varying  T.  For  small  T,  the 
period  of  the  stirring  is  less  than  the  timestep  of  integration  and  one  can  expect  the  flow  to 
appear  nondiffusive.  Indeed,  in  the  limit  of  small  T,  the  sine  map  given  in  (26)  approaches 
a  symplectic  integrator  [3].  This  result  can  be  obtained  by  noting  that 


x'  —  x 


T 


sin 


2irky\ 
— ) 


y  -y 


—  sm 


2itkx'  \ 

) 
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Figure  3:  Streamfunction  in  the  limit  that  T  — >•  0  .  The  flow  field  then  satisfies  u  = 
(dyip,  —dxip)  so  that  the  flow  is  clockwise  where  ip  is  negative  (dashed)  and  counterclockwise 
where  ip  is  positive  (solid). 


and  by  taking  the  limit  as  T  — >  0.  In  this  limit, 


^  etc.  such  that 


dx  (  2n rky\ 

—  =  sin  — - — 

dt  \  L  J 


dy  ( 2nkx 

—  =  sin  — - — 

dt  V  L 


(27) 


Equation  (27)  is  equivalent  to  saying  the  nondivergent  velocity  field,  u  =  dx/dt ,  can  be 
written  in  terms  of  a  streamfunction  u  =  (dyip,  —dxip)  where 


ip 


L 

2nk 


—  cos 


2irky\ 

L  ) 


(28) 


A  contour  plot  of  ip  is  given  in  figure  3. 

Consider  now  the  advection  of  a  cloud  of  particles  displaced  by  the  map  (26).  If  we  let 
L  =  k  =  1  and  use  a  relatively  small  value  of  T,  T  =  0.1,  the  flow  follows  nearly  closed 
orbits  as  given  in  figure  4.  As  this  flow  closely  approximates  the  streamfunction  solution, 
the  stirring  given  by  (26)  does  not  diffuse  the  concentration  of  particles.  As  a  result,  the 
use  of  an  ’effective’  diffusivity  does  not  apply  here. 

While  the  orbit  given  in  figure  4  is  not  strictly  confined  to  the  streamfunction  in  (28) 
because  of  the  use  of  a  finite  T,  it  is  not  space  filling.  We  can  add  noise  to  the  sinemap 
though 


x  =  x  +  T  sin 


+  \J2D  a 


y'  =  y  +  T  sin 


2tt  kx'\ 

) 


+  V2Db 


(29) 


where  a  and  b  are  Gaussian  random  variables  with  zero  mean  and  unit  variance.  Here  D  is  a 
constant  that  controls  the  amplitude  of  the  noise  and  can  be  consider  a  molecular  diffusivity. 
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itef  =  7,  T *  1  L=1  0  =  0 


(c)  D  =  0,  T  =  1,  and  it  =  7 


Figure  4:  Particle  concentration  in  a  box  of  sides  L  =  1  and  stirring  timestep  T  =  1  after 
20  iterations.  There  is  no  noise  in  (a)  and  the  map  is  confined  near  the  streamlines  given 
in  3.  With  the  addition  of  noise,  where  D  =  10-5,  the  orbit  becomes  diffusive.  Figure  4(c) 
is  the  same  as  in  4(a)  but  with  T  =  1  and  is  shown  only  at  the  7th  iteration. 


In  the  absence  of  stirring,  (29)  becomes  the  two-dimensional  version  of  the  random  walk 
model  given  in  (20).  As  we  have  chosen  (a2)  =  ( b 2)  =  1,  the  standard  deviation  a2  of  the 
noise  is  simply  2D. 

If  we  let  D  =  10~5  (figure  4(b)),  the  near-integrability  of  the  map  of  figure  4(a)  is 
broken,  however  the  flow  still  remains  constrained  to  near  the  streamlines  in  (28)  for  a 
moderate  number  of  iterations.  The  addition  of  this  molecular  diffusivity,  D ,  will  cause  the 
trajectory  to  eventually  fill  the  domain.  However,  this  process  is  quite  slow. 

The  diffusive  behavior  captured  by  adding  random  noise  in  figure  4(b)  can  also  be 
obtained  by  increasing  the  step  size  even  in  the  absence  of  noise.  For  a  box  size  L  =  1  and 
with  D  =  0  and  a  timestep  of  T  =  1,  ten  times  longer  than  that  of  the  examples  given 
in  figure  4,  the  flow  no  longer  follows  the  streamfunction  given  in  (28).  The  flow  in  figure 
4(c)  is  chaotic  and  quite  rapidly  fills  the  entire  domain  with  the  tracer  after  only  seven 
iterations. 

2.5  Filament  width  and  the  Batchelor  length 

The  Batchelor  length,  l  =  yj k/X ,  is  the  approximate  scale  at  which  filaments,  stretched 
by  stirring,  become  so  thin  that  the  tendency  for  the  filament  to  become  narrower  due  to 
stirring  is  balanced  by  the  diffusion  of  the  patch.  In  the  definition  of  the  Batchelor  length,  A 
is  the  local  rate  of  strain.  In  the  sine  map  given  in  (26),  modified  with  the  diffusion  equation 
of  (29),  the  rate  of  strain  is  approximately  T  and  the  coefficient  of  ”  molecular”  diffusion  is 
D.  For  this  problem,  then,  the  width  of  the  stirred  filaments  should  scale  approximately  as 

l  =  VlTT  (30) 

such  that  the  Batchelor  length  scales  with  \J~D  for  fixed  T . 

Snapshot  of  the  concentration  fields  for  ‘strong’  diffusivity,  with  D  =  10-4,  and  ‘weak’ 
diffusivity,  with  D  =  10~6,  are  given  in  5(a)  and  5(b),  respectively.  The  width  of  a  typical 
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0.5 


iter  =  3  T  =  1  L  =  1  D  :  0 0001 


iter  =  3.  T  =  1  L  =  1  D  =  1e-06 


Figure  5:  Particle  concentration  in  a  box  of  sides  L  =  1  with  T  =  1  and  with  varying  D  at 
iteration  3.  The  width  of  the  filaments  scale  as  the  Batchelor  length  given  in  (30). 


filament  in  the  ‘strong’  diffusivity  case  is  0.109  nondimensional  units  while  a  ‘weakly’  diffu¬ 
sive  filament  is  approximately  0.006  nondimensional  units.  The  width  of  the  filaments  for 
these  two  cases  do  roughly  scale  with  y/D  as  suggested  from  (30). 

2.6  Effective  diffusivity  and  the  domain  size,  L 

For  the  effective  diffusivity  scaling  given  in  (25)  to  hold  for  a  stirred  fluid,  the  length  scale 
of  the  stirring  must  be  less  than  that  of  the  domain.  In  the  examples  above,  the  length  scale 
of  the  stirring  was  the  same  as  that  of  the  domain.  If  we  increase  the  size  of  the  domain 
to  L  >  1  while  holding  T  at  an  intermediate  value,  T  =  1/2,  the  concept  of  an  ‘effective’ 
diffusivity  becomes  relevant. 

Shown  in  figure  6  are  snapshots  of  the  particle  concentration  field  after  ten  iterations  for 
a  domain  size  varying  from  L  =  1  to  L  =  25  while  the  spatial  correlations  ( x 2),  ( y 2),  and 
( X 2)  =  ( x 2  +  y 2)  and  cross-correlations  ( xy )  are  shown  in  figure  7.  For  small  domain  sizes, 
such  as  in  6(a)  and  6(b),  the  sine  flow  stirring  efficiently  fills  the  domain  with  the  tracer. 
This  homogenization  of  the  tracer  field  is  evident  in  6(c)  and  6(d)  by  the  leveling  off  of  the 
correlation  values  with  time. 

In  cases  with  a  relatively  small  domain  size,  such  as  6(a)  and  6(b),  the  spatial  variances 
for  ( x 2),  (y2,  and  ( X 2)  do  not  grow  linearly  with  time  as  predicted  from  the  effective 
diffusivity  equations  in  one  dimension,  equations  (22)  and  (23).  The  concept  of  an  effective 
diffusivity  fails  in  cases  where  the  stirring  is  at  approximately  the  same  scale  as  the  domain, 
as  predicted. 

Increasing  the  domain  size  to  L  =  10  and  L  =  25  while  holding  the  scale  of  the  stirring 
fixed,  however,  shows  a  more  promising  diffusive-like  behavior.  The  tracer  concentration  in 
6(c)  and  6(d)  spreads  outward  from  an  initial  patch  due  to  the  stirring  with  spatial  variances 
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Figure  6:  Particle  concentration  in  a  box  of  of  domain  size  of  L  =  1  in  6(a),  L  =  3  in  6(b), 
L  =  10  in  6(c),  and  L  =  25  in  6(d)  with  T  =  1  and  with  D=lx  10-6  at  iteration  10. 


growing  linearly  with  time,  as  shown  in  figures  7(c)  and  7(d). 

Note  that  the  variances  in  the  x  and  y  directions  for  all  domain  scales  indicate  isotropy,  as 
shown  by  the  similar  behavior  of  ( x 2)  and  ( y 2)  in  figure  6.  This  agrees  with  the  assumption 
that  the  variances  in  multiple  dimensions  have  the  same  statistics/standard  deviation  as  we 
assumed  in  the  derivation  of  the  multi-dimensional  form  of  the  effective  diffusivity  equation, 
(23).  Additionally,  the  cross-correlation  between  the  x  and  y  spatial  dimensions,  (xy),  is 
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Figure  7:  Calculations  of  the  spatial  variance  for  the  examples  shown  in  6  as  a  function  of 
iteration  and  increasing  domain  size. 


zero  at  long  times. 

From  (23),  we  calculate  the  approximate  effective  diffusivity  for  the  L  =  25  example 
as  At  =  0.068,  four  orders  of  magnitude  larger  than  the  molecular  diffusivity,  D  =  10~6  - 
the  stirring  in  this  case  more  efficiently  disperses  the  tracer  field  than  molecular  diffusivity 
alone. 


3  Matlab  Codes 

3.1  run.sinemap .m 


Listing  1:  run_sinemap.m 

function  [varargout]  =  run_ s inemap (T , L , D , Np , sqw) 

if  nargin  <  1,  T  =  1;  end 
if  nargin  <2,  L  =  1;  end 
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